Imports

library(evd); 
library(evdbayes);
library(coda);
library(ismev);
Lade nötiges Paket: mgcv
Lade nötiges Paket: nlme
This is mgcv 1.8-24. For overview type 'help("mgcv-package")'.
library(xts);
Lade nötiges Paket: zoo

Attache Paket: ‘zoo’

The following objects are masked from ‘package:base’:

    as.Date, as.Date.numeric

Loading the Data

load("../data/CAPE_Minder_Rychener_Malsot.RData")
load("../data/NINO34.RData")
load("../data/SRH_Minder_Rychener_Malsot.RData")

Generate PROD

prod = sqrt(cape)*srh

** Create Time Series Objects **

dates <- seq.Date(as.Date("1979-1-1"), as.Date("2015-12-31"), by=1)
feb29ix <- format(as.Date(dates), "%m-%d") == "02-29"
dates <- dates[!feb29ix]
prod_ts = xts(prod, order.by = rep(dates, each=8))

Beginning the Analysis

month_names = c("jan","feb","mar","apr","may","jun","jul","aug","sep","oct","nov","dec")
get_monthly = function(x) {
  output = list()
  len = nrow(x)
  
  # Get month of element
  month = time(x)
  month = gsub("....-", "", month)
  month = gsub("-..", "", month)
  monthlist = unique(month)
  for (i in 1:12) {
    output[[i]] = x[month == monthlist[i],]
  }
  names(output) = month_names
  return(output)
}
monthly_max = get_monthly(apply.monthly(prod_ts, max))
r = 2
r_monthly_max = get_monthly(apply.monthly(prod_ts, function(x) order(x, decreasing=T)[1:r]))
# get the monthly maxima of prod
m1 = as.data.frame(apply.monthly(prod_ts, max))$V1;
# produce the gumbel qq plot
gumbel_qq = function(x, title="") qqplot(qgumbel(c(1:length(x))/(length(x)+1)),
                                         x,
                                         main = paste("Gumbell Q-Q Plot", title),
                                         xlab = "Theoretical Quantiles", 
                                         ylab = "Sample Quantiles") ;
gumbel_qq(m1)

#qqplot(qgumbel(c(1:length(m1))/(length(m1)+1)),m1,main = "Gumbell Q-Q Plot",xlab = "Theoretical Quantiles", ylab = "Sample Quantiles") ;

The qq plot doesn’t fit very well, especially in the lower tail. This is likely due to seasonal dependence.

Fitting GEV to the entire data

# fit gevd with MLE and produce diagnostic plots
fitmax.MLE<-fgev(m1,method="Nelder-Mead")
par(mfrow=c(2,2))
fitmax.MLE

Call: fgev(x = m1, method = "Nelder-Mead") 
Deviance: 9272.599 

Estimates
      loc      scale      shape  
7.519e+03  7.013e+03  4.757e-03  

Standard Errors
      loc      scale      shape  
421.90201  331.43749    0.05347  

Optimization Information
  Convergence: successful 
  Function Evaluations: 171 
plot(fitmax.MLE)

Poor fit, probably because the distribution isn’t stationary. This is especially visible in the Probability plot, in which the confidence band is surpassed, indicating a poor fit.

# fit gevd with Bayesian Techniques
# use the previous outputs (rounded) as initial values (use a different shape)
init<-c(1.6e4,4e3,0.1)
# arbitrary priors
mat <- diag(c(10000,10000,100)) 
pn <- prior.norm(mean=c(0,0,0),cov=mat)
# proposal standard deviation (found by trial and error to get 20%<acceptance rate<40%)
psd<-c(500,0.03,0.02)
# draw 3k samples from markov chain
nit = 30000
thinning = 300
post<-posterior(nit, init=init,prior=pn,lh="gev",data=m1,psd=psd)
# diagnostic plots
MCMC<-mcmc(post[1:nit %% thinning == 0, ], thin=300) 
plot(MCMC) 

attr(mcmc(post),"ar")
            mu sigma   xi total
acc.rates 0.24  0.71 0.70  0.55
ext.rates 0.00  0.00 0.01  0.00
#MCMC_stab <- mcmc(post, thin=50, start=1000)
acf(mcmc(post[1:nit %% thinning == 0, ]))

We observe that there seem to be no substantial issues from the autocorrelation plots.

apply(mcmc(post[1:nit %% thinning == 0, ]),2,mean)
           mu         sigma            xi 
  236.0524535 11181.3806761    -0.1897774 
apply(mcmc(post[1:nit %% thinning == 0, ]),2,sd)
          mu        sigma           xi 
 98.60894002 624.46021143   0.03409059 

** Fit with MLE for months separately**

#monthly_fits = lapply(monthly_max, 
#                      function(x) fgev(data.frame(x)[,1], method="Nelder-Mead"))
monthly_fits = list()
error_cases = c(9, 12)
for (i in 1:length(monthly_max)) {
  print(i)
  
  if (i %in% error_cases) {
    monthly_fits[[i]] = fgev(as.data.frame(monthly_max[[i]])$V1, 
                             method = "Nelder-Mead",
                             std.err = FALSE)
  }
  else {
    monthly_fits[[i]] = fgev(as.data.frame(monthly_max[[i]])$V1, 
                             method = "Nelder-Mead")
  }
}
[1] 1
[1] 2
[1] 3
[1] 4
[1] 5
[1] 6
[1] 7
[1] 8
[1] 9
optimization may not have succeeded
[1] 10
[1] 11
[1] 12
names(monthly_fits) = names(monthly_max)
gumbel_qq(data.frame(monthly_max[[9]])[,1], "September")

gumbel_qq(data.frame(monthly_max[[12]])[,1], "December")

Fit R largest order statistics

largest_10 = get_monthly(apply.monthly(prod_ts, function(x) order(x, decreasing=T)[1:4]))
largest_2_fit = lapply(largest_10, function(x) rlarg.fit(x[ , 1:2]))
$conv
[1] 0

$nllh
[1] 381.7301

$mle
[1] 172.4434267  63.1813324  -0.8246329

$se
[1] 10.3518367  6.2427031  0.1133745

$conv
[1] 0

$nllh
[1] 366.2861

$mle
[1] 155.5451810  51.4464515  -0.7267874

$se
[1] 8.3567636 4.6071290 0.1019313

$conv
[1] 0

$nllh
[1] 348.4217

$mle
[1] 195.3919929  40.8441452  -0.8324715

$se
[1] 6.33288212 4.00545299 0.08820541

$conv
[1] 0

$nllh
[1] 375.033

$mle
[1] 162.9292046  57.8831308  -0.7755148

$se
[1] 9.282767 5.432561 0.099956

$conv
[1] 0

$nllh
[1] 380.4256

$mle
[1] 157.2651459  61.7442705  -0.6554637

$se
[1] 11.5888901  5.2748342  0.1508781

$conv
[1] 0

$nllh
[1] 364.4628

$mle
[1] 104.8006071  51.1390400  -0.1253138

$se
[1] 9.153534 4.541184 0.147210

$conv
[1] 0

$nllh
[1] 386.7889

$mle
[1] 128.6412577  66.3424033  -0.4883028

$se
[1] 11.877915  5.260236  0.139483

$conv
[1] 0

$nllh
[1] 380.7543

$mle
[1] 130.1181620  62.4551520  -0.3960438

$se
[1] 10.9040496  4.7896609  0.1237002

$conv
[1] 0

$nllh
[1] 361.8412

$mle
[1] 111.0849233  49.8496920  -0.2990619

$se
[1] 8.2812262 3.8324125 0.1013576

$conv
[1] 0

$nllh
[1] 377.358

$mle
[1] 75.1622227 56.9409902  0.2522703

$se
[1] 10.0789585  7.2986557  0.1984848

$conv
[1] 0

$nllh
[1] 376.1717

$mle
[1] 106.7510273  60.2707332  -0.1192325

$se
[1] 10.4205289  5.3759827  0.1333882

$conv
[1] 0

$nllh
[1] 377.2483

$mle
[1] 126.4665127  59.0621082  -0.2939312

$se
[1] 9.9836743 4.7019303 0.1170011
lapply(largest_2_fit, rlarg.diag)
Fehler in get(".xts_chob", .plotxtsEnv) : 
  Objekt '.xts_chob' nicht gefunden

get_se = function(x, ix) {
  if (is.null(x$std.err)) 0
  else x$std.err[ix]
}
mle_loc = unlist(lapply(monthly_fits, function(x) x$estimate[1]))
mle_loc_se = unlist(lapply(monthly_fits, get_se, 1))
mle_scale = unlist(lapply(monthly_fits, function(x) x$estimate[2]))
mle_scale_se = unlist(lapply(monthly_fits, get_se, 2))
mle_shape = unlist(lapply(monthly_fits, function(x) x$estimate[3]))
mle_shape_se = unlist(lapply(monthly_fits, get_se, 3))
plot_w_err = function(x, y, se, title = NULL) {
  max_ix = which.max(y)
  min_ix = which.min(y)
  plot(x, y,
       ylim = c(y[min_ix] - se[min_ix], y[max_ix] + se[max_ix]),
       main = title)
  arrows(x,y-se,x,y+se, code=3, length=0.02, angle = 90)
}
plot_w_err(1:12, mle_loc, mle_loc_se, "Location Parameter as Estimated with Likelihood")

plot_w_err(1:12, mle_scale, mle_scale_se, "Scale Parameter as Estimated with Likelihood")

plot_w_err(1:12, mle_shape, mle_shape_se, "Shape Parameter as Estimated with Likelihood")

** Fit with Bayesian Methods for months separately**

print(acceptance_rates)
$jan
   mu sigma    xi total 
 0.31  0.32  0.33  0.32 

$feb
   mu sigma    xi total 
 0.24  0.30  0.33  0.29 

$mar
   mu sigma    xi total 
 0.24  0.32  0.20  0.25 

$apr
   mu sigma    xi total 
 0.23  0.35  0.28  0.29 

$may
   mu sigma    xi total 
 0.24  0.29  0.21  0.25 

$jun
   mu sigma    xi total 
 0.24  0.32  0.22  0.26 

$jul
   mu sigma    xi total 
 0.23  0.20  0.29  0.24 

$aug
   mu sigma    xi total 
 0.24  0.21  0.25  0.23 

$sep
   mu sigma    xi total 
 0.23  0.27  0.20  0.23 

$oct
   mu sigma    xi total 
 0.24  0.38  0.23  0.28 

$nov
   mu sigma    xi total 
 0.24  0.34  0.34  0.31 

$dec
   mu sigma    xi total 
 0.34  0.33  0.30  0.33 

TODO-> R largest fit

PART 2 First, we check if the location parameter depends on time using a likelihood ratio test

#monthly_fits = lapply(monthly_max, 
#                      function(x) fgev(data.frame(x)[,1], method="Nelder-Mead"))
ratios = list()
trend = 1:length(as.data.frame(monthly_max[[12]])$V1)
trend = (trend-mean(trend))/sd(trend) # scale and center covariates as recommended
error_cases = c(9, 12)
for (i in 1:length(monthly_max)) {
  print(i)
  
  fit_const = fgev(as.data.frame(monthly_max[[i]])$V1, 
                             method = "Nelder-Mead",
                             std.err = FALSE)
  fit_dependant = fgev(as.data.frame(monthly_max[[i]])$V1, 
                             method = "Nelder-Mead",
                             nsloc = trend,
                             std.err = FALSE)
  
  ratios[[i]] = fit_const$dev-fit_dependant$dev 
}
[1] 1
optimization may not have succeeded
[1] 2
optimization may not have succeeded
[1] 3
[1] 4
optimization may not have succeeded
[1] 5
optimization may not have succeeded
[1] 6
[1] 7
[1] 8
[1] 9
optimization may not have succeeded
[1] 10
[1] 11
optimization may not have succeeded
[1] 12
names(ratios) = names(monthly_max)
chi_95level = qchisq(1-0.05/12,1)
plot(unlist(ratios),main="95% confidence test for time independance, Bonferroni multiple Testing", xlab="Month",ylab="Likelyhood Ratio Statistic")
abline(a=chi_95level,b=0,col="red")

Now, let’s check for independance from ENSO

#monthly_fits = lapply(monthly_max, 
#                      function(x) fgev(data.frame(x)[,1], method="Nelder-Mead"))
ratios = list()
# split nino data into months
n = nino34
dim(n)=c(12,length(as.data.frame(monthly_max[[12]])$V1))
error_cases = c(9, 12)
for (i in 1:length(monthly_max)) {
  print(i)
  trend = n[i,]
  trend = (trend-mean(trend))/sd(trend) # scale and center covariates as recommended
  fit_const = fgev(as.data.frame(monthly_max[[i]])$V1, 
                             method = "Nelder-Mead",
                             std.err = FALSE)
  fit_dependant = fgev(as.data.frame(monthly_max[[i]])$V1, 
                             method = "Nelder-Mead",
                             nsloc = trend,
                             std.err = FALSE)
  
  ratios[[i]] = fit_const$dev-fit_dependant$dev 
}
[1] 1
optimization may not have succeeded
[1] 2
optimization may not have succeeded
[1] 3
[1] 4
optimization may not have succeeded
[1] 5
[1] 6
optimization may not have succeeded
[1] 7
[1] 8
optimization may not have succeeded
[1] 9
optimization may not have succeeded
[1] 10
[1] 11
[1] 12
names(ratios) = names(monthly_max)
chi_95level = qchisq(1-0.05/12,1)
plot(unlist(ratios),main="95% confidence test for independance from ENSO, Bonferroni multiple testing", xlab="Month",ylab="Likelyhood Ratio Statistic")
abline(a=chi_95level,b=0,col="red")

Another method is the chi plot:

# plot the chi plot for dependance analysis
trend = 1:length(as.data.frame(monthly_max[[12]])$V1)
n = nino34
dim(n)=c(12,length(as.data.frame(monthly_max[[12]])$V1))
for (i in 1:length(monthly_max)) {
  nino = n[i,]
  m_data=as.data.frame(monthly_max[[i]])$V1
  dat.m1_month = cbind(m_data,trend);
  dat.m1_nino = cbind(m_data,nino);
  par(mfrow=c(2,2))
  chiplot(dat.m1_month,main1 = "Chi Plot Time",main2 = "Chi Bar Plot Time");
  chiplot(dat.m1_nino,main1 = "Chi Plot ENSO",main2 = "Chi Bar Plot ENSO");
}

PART 3 We will now analyse temporal clustering of extremes. For this, we will use the exiplot function from the evd library.

# define a function for getting the extremal indices for each month for a given threshold
monthly_eindex <- function(data, threshold_p, r=0){
  ei = list()
  for (i in 1:length(data)) {
    threshold = quantile(as.data.frame(data[[i]])$V1, threshold_p)
    ei[[i]]=exi(as.data.frame(data[[i]])$V1, threshold, r)
  }
  names(ei) = names(data)
  return(ei)
}
ei = monthly_eindex(get_monthly(prod_ts), 0.95)
plot(unlist(ei), main="Extremal Index by Month, 95%-Quantile as Threshold", xlab="Month", ylab="Extremal index")

We observe that the extremal index is 0.25-0.45, we can therefore conclude that we have strong dependance of extremes, with the limiting mean cluster size being roughly from 2 to 4. The clustering has no effect for estimaters based on the (monthly) maximum, but the r largest estimator is influenced by it.

PART 4 First, let’s estimate the 10 year return level using point process approach

monthly_fits_pp = list()
monthly_data = get_monthly(prod_ts)
error_cases = c(5,6,7,8,9,10,11,12)
month_days = c(31,28,31,30,31,30,31,31,30,31,30,31)
for (i in 1:length(monthly_max)) {
  print(i)
  threshold = quantile(as.data.frame(monthly_data[[i]])$V1, 0.95)
  
  if (i %in% error_cases) {
    monthly_fits_pp[[i]] = fpot(as.data.frame(monthly_data[[i]])$V1,
                             threshold = threshold,
                             model="pp",
                             npp = month_days[i]*8,
                             cmax = TRUE,
                             r = 1,
                             std.err = FALSE,
                             method = "Nelder-Mead")
  }
  else {
    monthly_fits_pp[[i]] = fpot(as.data.frame(monthly_data[[i]])$V1,
                             threshold = threshold,
                             model="pp",
                             npp = month_days[i]*8,
                             cmax = TRUE,
                             r = 1,
                             method = "Nelder-Mead")
  }
}
[1] 1
[1] 2
[1] 3
[1] 4
[1] 5
[1] 6
[1] 7
[1] 8
[1] 9
[1] 10
[1] 11
[1] 12
names(monthly_fits_pp) = names(monthly_data)
for(i in 1:12){
  par(mfrow=c(2,2)) 
  plot(monthly_fits_pp[[i]])
}

The fit in february has completely failed, and the others are not very good either

We will still estimate the return levels:

return_level = function(x,period=20){
  p = 1/period
  loc = x$estimate[[1]]
  scale = x$estimate[[2]]
  shape = x$estimate[[3]]
  level = loc + scale*(((-log(1-p))^-shape-1)/shape)
  return(level)
}
return_level_20 = lapply(monthly_fits, return_level) # 20 for testing
return_level_100 = lapply(monthly_fits, return_level, period=100)
return_level_50 = lapply(monthly_fits, return_level, period=50)
plot(unlist(return_level_100),main="100 Year Return level, estimated with point process", xlab="Month",ylab="Return Level")

plot(unlist(return_level_50),main="100 Year Return level, estimated with point process", xlab="Month",ylab="Return Level")

TODO: estimat with mcmc Assuming that we have the posterior densities of the markov chains, call theis function to plot return level histograms

return_level_mcmc = function(posterior,period=20){
  u = mc.quant(posterior,p=1-1/period,lh="gev")
  label_mcmc_rl = sprintf("%s Year Return Level",period)
  hist(u,nclass=20,prob=T,xlab=label_mcmc_rl, main = "Return Level Histogram")
}
lapply(monthly_bayes_fit, function(x) return_level_mcmc(x$posterior))

$jan
$breaks
 [1]     0  5000 10000 15000 20000 25000 30000 35000 40000 45000 50000 55000 60000 65000
[15] 70000 75000 80000 85000 90000

$counts
 [1] 14646 12954  1775   383   131    58    15    17     7     4     3     3     0     1
[15]     1     2     0     1

$density
 [1] 9.763675e-05 8.635712e-05 1.183294e-05 2.553248e-06 8.733042e-07 3.866538e-07
 [7] 9.999667e-08 1.133296e-07 4.666511e-08 2.666578e-08 1.999933e-08 1.999933e-08
[13] 0.000000e+00 6.666444e-09 6.666444e-09 1.333289e-08 0.000000e+00 6.666444e-09

$mids
 [1]  2500  7500 12500 17500 22500 27500 32500 37500 42500 47500 52500 57500 62500 67500
[15] 72500 77500 82500 87500

$xname
[1] "u"

$equidist
[1] TRUE

attr(,"class")
[1] "histogram"

$feb
$breaks
 [1]     0  5000 10000 15000 20000 25000 30000 35000 40000 45000 50000 55000 60000 65000

$counts
 [1]     1 18276 10221  1129   242    74    32    14     5     2     2     1     2

$density
 [1] 6.666444e-09 1.218359e-04 6.813773e-05 7.526416e-06 1.613280e-06 4.933169e-07
 [7] 2.133262e-07 9.333022e-08 3.333222e-08 1.333289e-08 1.333289e-08 6.666444e-09
[13] 1.333289e-08

$mids
 [1]  2500  7500 12500 17500 22500 27500 32500 37500 42500 47500 52500 57500 62500

$xname
[1] "u"

$equidist
[1] TRUE

attr(,"class")
[1] "histogram"

$mar
$breaks
 [1]  4000  6000  8000 10000 12000 14000 16000 18000 20000 22000 24000 26000 28000 30000
[15] 32000 34000 36000 38000 40000 42000 44000 46000 48000 50000 52000 54000

$counts
 [1]     1     1     0     0     0     0   103  1680  7830 11503  5998  1980   593   178
[15]    66    35    16    11     2     2     1     0     0     0     1

$density
 [1] 1.666611e-08 1.666611e-08 0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00
 [7] 1.716609e-06 2.799907e-05 1.304957e-04 1.917103e-04 9.996333e-05 3.299890e-05
[13] 9.883004e-06 2.966568e-06 1.099963e-06 5.833139e-07 2.666578e-07 1.833272e-07
[19] 3.333222e-08 3.333222e-08 1.666611e-08 0.000000e+00 0.000000e+00 0.000000e+00
[25] 1.666611e-08

$mids
 [1]  5000  7000  9000 11000 13000 15000 17000 19000 21000 23000 25000 27000 29000 31000
[15] 33000 35000 37000 39000 41000 43000 45000 47000 49000 51000 53000

$xname
[1] "u"

$equidist
[1] TRUE

attr(,"class")
[1] "histogram"

$apr
$breaks
 [1]     0  5000 10000 15000 20000 25000 30000 35000 40000 45000 50000 55000 60000 65000
[15] 70000

$counts
 [1]     1     3     3     3    88 20707  8892   270    21     7     1     3     1     1

$density
 [1] 6.666444e-09 1.999933e-08 1.999933e-08 1.999933e-08 5.866471e-07 1.380421e-04
 [7] 5.927802e-05 1.799940e-06 1.399953e-07 4.666511e-08 6.666444e-09 1.999933e-08
[13] 6.666444e-09 6.666444e-09

$mids
 [1]  2500  7500 12500 17500 22500 27500 32500 37500 42500 47500 52500 57500 62500 67500

$xname
[1] "u"

$equidist
[1] TRUE

attr(,"class")
[1] "histogram"

$may
$breaks
 [1]     0  5000 10000 15000 20000 25000 30000 35000 40000 45000 50000 55000 60000 65000
[15] 70000 75000 80000 85000

$counts
 [1]     1     3     3     3     1     1  1979 24216  3640   125    15     7     2     2
[15]     2     0     1

$density
 [1] 6.666444e-09 1.999933e-08 1.999933e-08 1.999933e-08 6.666444e-09 6.666444e-09
 [7] 1.319289e-05 1.614346e-04 2.426586e-05 8.333056e-07 9.999667e-08 4.666511e-08
[13] 1.333289e-08 1.333289e-08 1.333289e-08 0.000000e+00 6.666444e-09

$mids
 [1]  2500  7500 12500 17500 22500 27500 32500 37500 42500 47500 52500 57500 62500 67500
[15] 72500 77500 82500

$xname
[1] "u"

$equidist
[1] TRUE

attr(,"class")
[1] "histogram"

$jun
$breaks
 [1]     0  5000 10000 15000 20000 25000 30000 35000 40000 45000 50000 55000 60000 65000
[15] 70000 75000 80000 85000

$counts
 [1]     1     3     3     3     1     1  2192 22381  5201   185    17     6     2     2
[15]     2     0     1

$density
 [1] 6.666444e-09 1.999933e-08 1.999933e-08 1.999933e-08 6.666444e-09 6.666444e-09
 [7] 1.461285e-05 1.492017e-04 3.467218e-05 1.233292e-06 1.133296e-07 3.999867e-08
[13] 1.333289e-08 1.333289e-08 1.333289e-08 0.000000e+00 6.666444e-09

$mids
 [1]  2500  7500 12500 17500 22500 27500 32500 37500 42500 47500 52500 57500 62500 67500
[15] 72500 77500 82500

$xname
[1] "u"

$equidist
[1] TRUE

attr(,"class")
[1] "histogram"

$jul
$breaks
 [1]  4000  6000  8000 10000 12000 14000 16000 18000 20000 22000 24000 26000 28000 30000
[15] 32000 34000 36000 38000 40000 42000

$counts
 [1]     2     3     1     3     2     0     1     0   522 13995 14730   663    54    15
[15]     6     2     1     0     1

$density
 [1] 3.333222e-08 4.999833e-08 1.666611e-08 4.999833e-08 3.333222e-08 0.000000e+00
 [7] 1.666611e-08 0.000000e+00 8.699710e-06 2.332422e-04 2.454918e-04 1.104963e-05
[13] 8.999700e-07 2.499917e-07 9.999667e-08 3.333222e-08 1.666611e-08 0.000000e+00
[19] 1.666611e-08

$mids
 [1]  5000  7000  9000 11000 13000 15000 17000 19000 21000 23000 25000 27000 29000 31000
[15] 33000 35000 37000 39000 41000

$xname
[1] "u"

$equidist
[1] TRUE

attr(,"class")
[1] "histogram"

$aug
$breaks
 [1]  4000  6000  8000 10000 12000 14000 16000 18000 20000 22000 24000 26000 28000 30000
[15] 32000 34000 36000 38000 40000 42000 44000

$counts
 [1]     1     1     0     0     2     0     0     6   216  3344 12190 10653  2914   518
[15]   115    20    14     2     2     3

$density
 [1] 1.666611e-08 1.666611e-08 0.000000e+00 0.000000e+00 3.333222e-08 0.000000e+00
 [7] 0.000000e+00 9.999667e-08 3.599880e-06 5.573148e-05 2.031599e-04 1.775441e-04
[13] 4.856505e-05 8.633046e-06 1.916603e-06 3.333222e-07 2.333256e-07 3.333222e-08
[19] 3.333222e-08 4.999833e-08

$mids
 [1]  5000  7000  9000 11000 13000 15000 17000 19000 21000 23000 25000 27000 29000 31000
[15] 33000 35000 37000 39000 41000 43000

$xname
[1] "u"

$equidist
[1] TRUE

attr(,"class")
[1] "histogram"

$sep
$breaks
 [1]     0  5000 10000 15000 20000 25000 30000 35000 40000 45000 50000 55000 60000 65000

$counts
 [1]     1     3     3    15 23255  6680    32     4     0     2     2     3     1

$density
 [1] 6.666444e-09 1.999933e-08 1.999933e-08 9.999667e-08 1.550282e-04 4.453185e-05
 [7] 2.133262e-07 2.666578e-08 0.000000e+00 1.333289e-08 1.333289e-08 1.999933e-08
[13] 6.666444e-09

$mids
 [1]  2500  7500 12500 17500 22500 27500 32500 37500 42500 47500 52500 57500 62500

$xname
[1] "u"

$equidist
[1] TRUE

attr(,"class")
[1] "histogram"

$oct
$breaks
 [1]  4000  6000  8000 10000 12000 14000 16000 18000 20000 22000 24000 26000 28000 30000
[15] 32000 34000 36000 38000 40000 42000 44000 46000 48000

$counts
 [1]    1    1    0    0    0   30  938 5430 9239 8150 3932 1437  485  204   70   44   20
[18]    9    3    3    1    4

$density
 [1] 1.666611e-08 1.666611e-08 0.000000e+00 0.000000e+00 0.000000e+00 4.999833e-07
 [7] 1.563281e-05 9.049698e-05 1.539782e-04 1.358288e-04 6.553115e-05 2.394920e-05
[13] 8.083064e-06 3.399887e-06 1.166628e-06 7.333089e-07 3.333222e-07 1.499950e-07
[19] 4.999833e-08 4.999833e-08 1.666611e-08 6.666444e-08

$mids
 [1]  5000  7000  9000 11000 13000 15000 17000 19000 21000 23000 25000 27000 29000 31000
[15] 33000 35000 37000 39000 41000 43000 45000 47000

$xname
[1] "u"

$equidist
[1] TRUE

attr(,"class")
[1] "histogram"

$nov
$breaks
 [1]  4000  6000  8000 10000 12000 14000 16000 18000 20000 22000 24000 26000 28000 30000
[15] 32000 34000 36000 38000 40000 42000 44000

$counts
 [1]     1   102  5217 14247  7181  2090   689   269   101    50    19    14     8     1
[15]     2     5     2     2     0     1

$density
 [1] 1.666611e-08 1.699943e-06 8.694710e-05 2.374421e-04 1.196793e-04 3.483217e-05
 [7] 1.148295e-05 4.483184e-06 1.683277e-06 8.333056e-07 3.166561e-07 2.333256e-07
[13] 1.333289e-07 1.666611e-08 3.333222e-08 8.333056e-08 3.333222e-08 3.333222e-08
[19] 0.000000e+00 1.666611e-08

$mids
 [1]  5000  7000  9000 11000 13000 15000 17000 19000 21000 23000 25000 27000 29000 31000
[15] 33000 35000 37000 39000 41000 43000

$xname
[1] "u"

$equidist
[1] TRUE

attr(,"class")
[1] "histogram"

$dec
$breaks
 [1]     0  5000 10000 15000 20000 25000 30000 35000 40000 45000 50000 55000 60000 65000
[15] 70000

$counts
 [1] 23735  5642   470   103    16    18     9     1     2     0     0     0     3     2

$density
 [1] 1.582281e-04 3.761208e-05 3.133229e-06 6.866438e-07 1.066631e-07 1.199960e-07
 [7] 5.999800e-08 6.666444e-09 1.333289e-08 0.000000e+00 0.000000e+00 0.000000e+00
[13] 1.999933e-08 1.333289e-08

$mids
 [1]  2500  7500 12500 17500 22500 27500 32500 37500 42500 47500 52500 57500 62500 67500

$xname
[1] "u"

$equidist
[1] TRUE

attr(,"class")
[1] "histogram"

LS0tCnRpdGxlOiAiVGh1bmRlcnN0cm9tIEFuYWx5c2lzIgpvdXRwdXQ6IGh0bWxfbm90ZWJvb2sKLS0tCioqSW1wb3J0cyoqCmBgYHtyfQpsaWJyYXJ5KGV2ZCk7IApsaWJyYXJ5KGV2ZGJheWVzKTsKbGlicmFyeShjb2RhKTsKbGlicmFyeShpc21ldik7CmxpYnJhcnkoeHRzKTsKCmBgYAoKCioqTG9hZGluZyB0aGUgRGF0YSoqCgpgYGB7cn0KbG9hZCgiLi4vZGF0YS9DQVBFX01pbmRlcl9SeWNoZW5lcl9NYWxzb3QuUkRhdGEiKQpsb2FkKCIuLi9kYXRhL05JTk8zNC5SRGF0YSIpCmxvYWQoIi4uL2RhdGEvU1JIX01pbmRlcl9SeWNoZW5lcl9NYWxzb3QuUkRhdGEiKQpgYGAKCioqR2VuZXJhdGUgUFJPRCoqCmBgYHtyfQpwcm9kID0gc3FydChjYXBlKSpzcmgKYGBgCgoKCgoqKiBDcmVhdGUgVGltZSBTZXJpZXMgT2JqZWN0cyAqKgpgYGB7cn0KZGF0ZXMgPC0gc2VxLkRhdGUoYXMuRGF0ZSgiMTk3OS0xLTEiKSwgYXMuRGF0ZSgiMjAxNS0xMi0zMSIpLCBieT0xKQpmZWIyOWl4IDwtIGZvcm1hdChhcy5EYXRlKGRhdGVzKSwgIiVtLSVkIikgPT0gIjAyLTI5IgpkYXRlcyA8LSBkYXRlc1shZmViMjlpeF0KCnByb2RfdHMgPSB4dHMocHJvZCwgb3JkZXIuYnkgPSByZXAoZGF0ZXMsIGVhY2g9OCkpCmBgYAoKCgoqKkJlZ2lubmluZyB0aGUgQW5hbHlzaXMqKgpgYGB7cn0KbW9udGhfbmFtZXMgPSBjKCJqYW4iLCJmZWIiLCJtYXIiLCJhcHIiLCJtYXkiLCJqdW4iLCJqdWwiLCJhdWciLCJzZXAiLCJvY3QiLCJub3YiLCJkZWMiKQpnZXRfbW9udGhseSA9IGZ1bmN0aW9uKHgpIHsKICBvdXRwdXQgPSBsaXN0KCkKICBsZW4gPSBucm93KHgpCiAgCiAgIyBHZXQgbW9udGggb2YgZWxlbWVudAogIG1vbnRoID0gdGltZSh4KQogIG1vbnRoID0gZ3N1YigiLi4uLi0iLCAiIiwgbW9udGgpCiAgbW9udGggPSBnc3ViKCItLi4iLCAiIiwgbW9udGgpCiAgbW9udGhsaXN0ID0gdW5pcXVlKG1vbnRoKQogIGZvciAoaSBpbiAxOjEyKSB7CiAgICBvdXRwdXRbW2ldXSA9IHhbbW9udGggPT0gbW9udGhsaXN0W2ldLF0KICB9CiAgbmFtZXMob3V0cHV0KSA9IG1vbnRoX25hbWVzCiAgcmV0dXJuKG91dHB1dCkKfQoKbW9udGhseV9tYXggPSBnZXRfbW9udGhseShhcHBseS5tb250aGx5KHByb2RfdHMsIG1heCkpCnIgPSAyCnJfbW9udGhseV9tYXggPSBnZXRfbW9udGhseShhcHBseS5tb250aGx5KHByb2RfdHMsIGZ1bmN0aW9uKHgpIG9yZGVyKHgsIGRlY3JlYXNpbmc9VClbMTpyXSkpCmBgYAoKCmBgYHtyfQojIGdldCB0aGUgbW9udGhseSBtYXhpbWEgb2YgcHJvZAptMSA9IGFzLmRhdGEuZnJhbWUoYXBwbHkubW9udGhseShwcm9kX3RzLCBtYXgpKSRWMTsKIyBwcm9kdWNlIHRoZSBndW1iZWwgcXEgcGxvdApndW1iZWxfcXEgPSBmdW5jdGlvbih4LCB0aXRsZT0iIikgcXFwbG90KHFndW1iZWwoYygxOmxlbmd0aCh4KSkvKGxlbmd0aCh4KSsxKSksCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgeCwKICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICBtYWluID0gcGFzdGUoIkd1bWJlbGwgUS1RIFBsb3QiLCB0aXRsZSksCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgeGxhYiA9ICJUaGVvcmV0aWNhbCBRdWFudGlsZXMiLCAKICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICB5bGFiID0gIlNhbXBsZSBRdWFudGlsZXMiKSA7CgpndW1iZWxfcXEobTEpCgojcXFwbG90KHFndW1iZWwoYygxOmxlbmd0aChtMSkpLyhsZW5ndGgobTEpKzEpKSxtMSxtYWluID0gIkd1bWJlbGwgUS1RIFBsb3QiLHhsYWIgPSAiVGhlb3JldGljYWwgUXVhbnRpbGVzIiwgeWxhYiA9ICJTYW1wbGUgUXVhbnRpbGVzIikgOwpgYGAKVGhlIHFxIHBsb3QgZG9lc24ndCBmaXQgdmVyeSB3ZWxsLCBlc3BlY2lhbGx5IGluIHRoZSBsb3dlciB0YWlsLiBUaGlzIGlzIGxpa2VseSBkdWUKdG8gc2Vhc29uYWwgZGVwZW5kZW5jZS4KCioqRml0dGluZyBHRVYgdG8gdGhlIGVudGlyZSBkYXRhKioKYGBge3J9CiMgZml0IGdldmQgd2l0aCBNTEUgYW5kIHByb2R1Y2UgZGlhZ25vc3RpYyBwbG90cwpmaXRtYXguTUxFPC1mZ2V2KG0xLG1ldGhvZD0iTmVsZGVyLU1lYWQiKQpwYXIobWZyb3c9YygyLDIpKQpmaXRtYXguTUxFCnBsb3QoZml0bWF4Lk1MRSkKYGBgClBvb3IgZml0LCBwcm9iYWJseSBiZWNhdXNlIHRoZSBkaXN0cmlidXRpb24gaXNuJ3Qgc3RhdGlvbmFyeS4gVGhpcyBpcyBlc3BlY2lhbGx5IAp2aXNpYmxlIGluIHRoZSBQcm9iYWJpbGl0eSBwbG90LCBpbiB3aGljaCB0aGUgY29uZmlkZW5jZSBiYW5kIGlzIHN1cnBhc3NlZCwgCmluZGljYXRpbmcgYSBwb29yIGZpdC4KCgpgYGB7cn0KIyBmaXQgZ2V2ZCB3aXRoIEJheWVzaWFuIFRlY2huaXF1ZXMKIyB1c2UgdGhlIHByZXZpb3VzIG91dHB1dHMgKHJvdW5kZWQpIGFzIGluaXRpYWwgdmFsdWVzICh1c2UgYSBkaWZmZXJlbnQgc2hhcGUpCmluaXQ8LWMoMS42ZTQsNGUzLDAuMSkKIyBhcmJpdHJhcnkgcHJpb3JzCm1hdCA8LSBkaWFnKGMoMTAwMDAsMTAwMDAsMTAwKSkgCnBuIDwtIHByaW9yLm5vcm0obWVhbj1jKDAsMCwwKSxjb3Y9bWF0KQojIHByb3Bvc2FsIHN0YW5kYXJkIGRldmlhdGlvbiAoZm91bmQgYnkgdHJpYWwgYW5kIGVycm9yIHRvIGdldCAyMCU8YWNjZXB0YW5jZSByYXRlPDQwJSkKcHNkPC1jKDUwMCwwLjAzLDAuMDIpCiMgZHJhdyAzayBzYW1wbGVzIGZyb20gbWFya292IGNoYWluCm5pdCA9IDMwMDAwCnRoaW5uaW5nID0gMzAwCnBvc3Q8LXBvc3RlcmlvcihuaXQsIGluaXQ9aW5pdCxwcmlvcj1wbixsaD0iZ2V2IixkYXRhPW0xLHBzZD1wc2QpCiMgZGlhZ25vc3RpYyBwbG90cwpNQ01DPC1tY21jKHBvc3RbMTpuaXQgJSUgdGhpbm5pbmcgPT0gMCwgXSwgdGhpbj0zMDApIApwbG90KE1DTUMpIAphdHRyKG1jbWMocG9zdCksImFyIikKCmBgYAoKCmBgYHtyfQojTUNNQ19zdGFiIDwtIG1jbWMocG9zdCwgdGhpbj01MCwgc3RhcnQ9MTAwMCkKYWNmKG1jbWMocG9zdFsxOm5pdCAlJSB0aGlubmluZyA9PSAwLCBdKSkKYGBgCldlIG9ic2VydmUgdGhhdCB0aGVyZSBzZWVtIHRvIGJlIG5vIHN1YnN0YW50aWFsIGlzc3VlcyBmcm9tIHRoZSBhdXRvY29ycmVsYXRpb24gCnBsb3RzLiAKCmBgYHtyfQphcHBseShtY21jKHBvc3RbMTpuaXQgJSUgdGhpbm5pbmcgPT0gMCwgXSksMixtZWFuKQphcHBseShtY21jKHBvc3RbMTpuaXQgJSUgdGhpbm5pbmcgPT0gMCwgXSksMixzZCkKCmBgYAoKKiogRml0IHdpdGggTUxFIGZvciBtb250aHMgc2VwYXJhdGVseSoqCmBgYHtyfQojbW9udGhseV9maXRzID0gbGFwcGx5KG1vbnRobHlfbWF4LCAKIyAgICAgICAgICAgICAgICAgICAgICBmdW5jdGlvbih4KSBmZ2V2KGRhdGEuZnJhbWUoeClbLDFdLCBtZXRob2Q9Ik5lbGRlci1NZWFkIikpCm1vbnRobHlfZml0cyA9IGxpc3QoKQplcnJvcl9jYXNlcyA9IGMoOSwgMTIpCmZvciAoaSBpbiAxOmxlbmd0aChtb250aGx5X21heCkpIHsKICBwcmludChpKQogIAogIGlmIChpICVpbiUgZXJyb3JfY2FzZXMpIHsKICAgIG1vbnRobHlfZml0c1tbaV1dID0gZmdldihhcy5kYXRhLmZyYW1lKG1vbnRobHlfbWF4W1tpXV0pJFYxLCAKICAgICAgICAgICAgICAgICAgICAgICAgICAgICBtZXRob2QgPSAiTmVsZGVyLU1lYWQiLAogICAgICAgICAgICAgICAgICAgICAgICAgICAgIHN0ZC5lcnIgPSBGQUxTRSkKICB9CiAgZWxzZSB7CiAgICBtb250aGx5X2ZpdHNbW2ldXSA9IGZnZXYoYXMuZGF0YS5mcmFtZShtb250aGx5X21heFtbaV1dKSRWMSwgCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgbWV0aG9kID0gIk5lbGRlci1NZWFkIikKICB9Cn0KbmFtZXMobW9udGhseV9maXRzKSA9IG5hbWVzKG1vbnRobHlfbWF4KQoKYGBgCgpgYGB7cn0KZ3VtYmVsX3FxKGRhdGEuZnJhbWUobW9udGhseV9tYXhbWzldXSlbLDFdLCAiU2VwdGVtYmVyIikKZ3VtYmVsX3FxKGRhdGEuZnJhbWUobW9udGhseV9tYXhbWzEyXV0pWywxXSwgIkRlY2VtYmVyIikKYGBgCgoqKkZpdCBSIGxhcmdlc3Qgb3JkZXIgc3RhdGlzdGljcyoqCmBgYHtyfQoKbGFyZ2VzdF8xMCA9IGdldF9tb250aGx5KGFwcGx5Lm1vbnRobHkocHJvZF90cywgZnVuY3Rpb24oeCkgb3JkZXIoeCwgZGVjcmVhc2luZz1UKVsxOjRdKSkKbGFyZ2VzdF8yX2ZpdCA9IGxhcHBseShsYXJnZXN0XzEwLCBmdW5jdGlvbih4KSBybGFyZy5maXQoeFsgLCAxOjJdKSkKYGBgCgpgYGB7cn0KZ2V0X3NlID0gZnVuY3Rpb24oeCwgaXgpIHsKICBpZiAoaXMubnVsbCh4JHN0ZC5lcnIpKSAwCiAgZWxzZSB4JHN0ZC5lcnJbaXhdCn0KbWxlX2xvYyA9IHVubGlzdChsYXBwbHkobW9udGhseV9maXRzLCBmdW5jdGlvbih4KSB4JGVzdGltYXRlWzFdKSkKbWxlX2xvY19zZSA9IHVubGlzdChsYXBwbHkobW9udGhseV9maXRzLCBnZXRfc2UsIDEpKQptbGVfc2NhbGUgPSB1bmxpc3QobGFwcGx5KG1vbnRobHlfZml0cywgZnVuY3Rpb24oeCkgeCRlc3RpbWF0ZVsyXSkpCm1sZV9zY2FsZV9zZSA9IHVubGlzdChsYXBwbHkobW9udGhseV9maXRzLCBnZXRfc2UsIDIpKQptbGVfc2hhcGUgPSB1bmxpc3QobGFwcGx5KG1vbnRobHlfZml0cywgZnVuY3Rpb24oeCkgeCRlc3RpbWF0ZVszXSkpCm1sZV9zaGFwZV9zZSA9IHVubGlzdChsYXBwbHkobW9udGhseV9maXRzLCBnZXRfc2UsIDMpKQpgYGAKCgpgYGB7cn0KcGxvdF93X2VyciA9IGZ1bmN0aW9uKHgsIHksIHNlLCB0aXRsZSA9IE5VTEwpIHsKICBtYXhfaXggPSB3aGljaC5tYXgoeSkKICBtaW5faXggPSB3aGljaC5taW4oeSkKICBwbG90KHgsIHksCiAgICAgICB5bGltID0gYyh5W21pbl9peF0gLSBzZVttaW5faXhdLCB5W21heF9peF0gKyBzZVttYXhfaXhdKSwKICAgICAgIG1haW4gPSB0aXRsZSkKICBhcnJvd3MoeCx5LXNlLHgseStzZSwgY29kZT0zLCBsZW5ndGg9MC4wMiwgYW5nbGUgPSA5MCkKfQpwbG90X3dfZXJyKDE6MTIsIG1sZV9sb2MsIG1sZV9sb2Nfc2UsICJMb2NhdGlvbiBQYXJhbWV0ZXIgYXMgRXN0aW1hdGVkIHdpdGggTGlrZWxpaG9vZCIpCnBsb3Rfd19lcnIoMToxMiwgbWxlX3NjYWxlLCBtbGVfc2NhbGVfc2UsICJTY2FsZSBQYXJhbWV0ZXIgYXMgRXN0aW1hdGVkIHdpdGggTGlrZWxpaG9vZCIpCnBsb3Rfd19lcnIoMToxMiwgbWxlX3NoYXBlLCBtbGVfc2hhcGVfc2UsICJTaGFwZSBQYXJhbWV0ZXIgYXMgRXN0aW1hdGVkIHdpdGggTGlrZWxpaG9vZCIpCgpgYGAKCioqIEZpdCB3aXRoIEJheWVzaWFuIE1ldGhvZHMgZm9yIG1vbnRocyBzZXBhcmF0ZWx5KioKYGBge3J9CgoKIyBGaXRzIEdFViBkaXN0cmlidXRpb24gd2l0aCBiYXllc2lhbiBtZXRob2QgZm9yIGdpdmVuIHBhcmFtZXRlcnMKYmF5ZXNfZml0dGVyID0gZnVuY3Rpb24oeCwgCiAgICAgICAgICAgICAgICAgICAgICAgIGluaXQgPSBjKDFlMywgMWUzLCAwLjEpLCAjIEluaXRpYWwgdmFsdWVzCiAgICAgICAgICAgICAgICAgICAgICAgIG1hdCA9IGRpYWcoYygxMDAwMCwxMDAwMCwxMDApKSwKICAgICAgICAgICAgICAgICAgICAgICAgcHNkID0gYyg1MDAsMC4xLDAuMSksICMgUHJvcG9zZWQgU0RldgogICAgICAgICAgICAgICAgICAgICAgICBuaXQgPSAzMDAwLCAjIE5iIEl0ZXJhdGlvbnMKICAgICAgICAgICAgICAgICAgICAgICAgdGhpbm5pbmcgPSA1MCwgIyBUaGlubmluZyBGYWN0b3IKICAgICAgICAgICAgICAgICAgICAgICAgZG9fZGlhZ24gPSBGQUxTRSwgIyBCb29sIHdoZXRoZXIgdG8gc2hvdyBkaWFnbm9zdGljIHBsb3RzCiAgICAgICAgICAgICAgICAgICAgICAgIGRvX2F1dG9yZWcgPSBGQUxTRSwgIyBCb29sIHdoZXRoZXIgdG8gc2hvdyBhdXRvcmVnIHBsb3RzCiAgICAgICAgICAgICAgICAgICAgICAgIHNlZWQgPSA0MiAjIFNlZWQgCiAgICAgICAgICAgICAgICAgICAgICAgICkgewogIHNldC5zZWVkKHNlZWQpICAgICAgICAgICAgICAgIAogIHBuID0gcHJpb3Iubm9ybShtZWFuPWMoMCwwLDApLGNvdj1tYXQpCiAgcG9zdCA9IHBvc3RlcmlvcihuaXQsIGluaXQ9aW5pdCwgcHJpb3I9cG4sIGxoPSJnZXYiLCBkYXRhPXgsIHBzZD1wc2QpCiAgCiAgaWYoZG9fZGlhZ24pIHsKICAgIE1DTUM8LW1jbWMocG9zdFsxOm5pdCAlJSB0aGlubmluZyA9PSAwLCBdKSAKICAgIHBsb3QoTUNNQykgCiAgfQogIGlmKGRvX2F1dG9yZWcpIHsKICAgIGFjZihtY21jKHBvc3RbMTpuaXQgJSUgdGhpbm5pbmcgPT0gMCwgXSkpCiAgfQogIGxpc3QocG9zdGVyaW9yID0gcG9zdCwgCiAgICAgICBhY2NlcHRhbmNlX3JhdGUgPSBhdHRyKG1jbWMocG9zdCksImFyIikpCn0KCgoKIyBJdGVyYXRpdmVseSBmaXRzIEdFViB3aXRoIGJheWVzaWFuIG1ldGhvZHMsIHVudGlsIHRoZSBmaXQgaGFzIAojIGFjY2VwdGFibGUgYWNjZXB0YW5jZSByYXRlcyAoaS5lLiAwLjIgPCBBUiA8IDAuNCkuIElmIHRoZSBBUiBpcyB0b28gaGlnaCwgCiMgdGhlIHByb3Bvc2VkIFNEZXYgZm9yIHRoZSBwYXJhbWV0ZXIgaXMgbXVsdGlwbGllZCB3aXRoIDEuNS4gSWYgaXQncyB0b28gc21hbGwsCiMgdGhlIHByb3Bvc2VkIFNEZXYgaXMgZGl2aWRlZCBieSAyLiBUaGlzIGlzIHJlcGVhdGVkIHVudGlsIHRoZSBhY2NlcHRhbmNlIHJhdGUKIyBpcyBnb29kIGZvciBhbGwgcGFyYW1ldGVycywgb3IgbWF4X2l0IGlzIHJlYWNoZWQuIFRoZW4sIGEgZmluYWwgbW9kZWwgaXMgZml0dGVkCiMgd2l0aCBtb3JlIGl0ZXJhdGlvbnMuIApiYXllc19maXR0ZXJfcGFyYW1fc2VhcmNoID0gZnVuY3Rpb24oeCwKICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIHBzZF9pbml0ID0gYyg1MDAsMC4xLDAuMSksICMgSW5pdGlhbCBwcm9wb3NlZCBTRGV2CiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICBuaXRfZnVsbCA9IDMwMDAsICMgTmIgSXRlcmF0aW9ucyBmb3IgZmluYWwgbW9kZWwKICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIG5pdF9zZWFyY2ggPSAxNTAsICMgTmIgSXRlcmF0aW9ucyBmb3IgcGFyYW0gc2VhcmNoCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICBkb19kaWFnbiA9IEZBTFNFLCAjIEJvb2wgd2hldGhlciB0byBzaG93IGRpYWdub3N0aWMgcGxvdHMKICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIGRvX2F1dG9yZWcgPSBGQUxTRSwgIyBCb29sIHdoZXRoZXIgdG8gc2hvdyBhdXRvcmVnIHBsb3RzCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICBtYXhfaXQgPSAyMCwKICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIC4uLiAjIEFkZGl0aW9uYWwgcGFyYW1zIHBhc3NlZCB0byBiYXllc19maXR0ZXIKICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICkgCnsKICAjIEl0ZXJhdGUgdW50aWwgZGVzaXJlZCBhY2NlcHRhbmNlIHJhdGUgaXMgb2J0YWluZWQKICBjb250ID0gVFJVRQogIHBzZCA9IHBzZF9pbml0CiAgaXQgPSAwCiAgd2hpbGUoY29udCkgewogICAgaXQgPSBpdCsxCiAgICBpZiAoaXQgPiBtYXhfaXQpIHsKICAgICAgd2FybmluZygiVGhlICIpCiAgICB9CiAgICBmaXQgPSBiYXllc19maXR0ZXIoeCwgcHNkPXBzZCwgbml0PW5pdF9zZWFyY2gsIGRvX2RpYWduPUZBTFNFLCAKICAgICAgICAgICAgICAgICAgICAgICBkb19hdXRvcmVnPUZBTFNFLC4uLikKICAgIGFjY19yYXRlcyA9IGZpdCRhY2NlcHRhbmNlX3JhdGVbMSwgMTozXQogICAgCiAgICB0b29faGlnaCA9IGFjY19yYXRlcyA+IC40CiAgICB0b29fbG93ID0gYWNjX3JhdGVzIDwgLjIKICAgIAogICAgaWYgKGFsbCghdG9vX2hpZ2gpICYmIGFsbCghdG9vX2xvdykpIHsgIyBBbGwgYWNjZXB0YW5jZSByYXRlcyBsaWUgd2l0aGluIHRocmVzaG9sZAogICAgICBjb250PUZBTFNFCiAgICB9IGVsc2UgaWYgKGl0ID4gbWF4X2l0KSB7ICMgbWF4X2l0IGlzIHJlYWNoZWQKICAgICAgY29udD1GQUxTRQogICAgICB3YXJuaW5nKCJtYXhfaXQgd2FzIHJlYWNoZWQiKQogICAgfSBlbHNlIHsgIyBDb3JyZWN0IHZhbHVlcyB3aGljaCBoYXZlIHdyb25nIHRocmVzaG9sZAogICAgICBwc2RbdG9vX2hpZ2hdID0gcHNkW3Rvb19oaWdoXSAqIDEuNQogICAgICBwc2RbdG9vX2xvd10gPSBwc2RbdG9vX2xvd10gLyAyCiAgICB9CiAgfQogIAogICMgRml0IGZpbmFsIG1vZGVsCiAgYmF5ZXNfZml0dGVyKHgsIHBzZD1wc2QsIG5pdD1uaXRfZnVsbCwgZG9fZGlhZ249ZG9fZGlhZ24sIAogICAgICAgICAgICAgICBkb19hdXRvcmVnPWRvX2F1dG9yZWcsIC4uLikKICAKfQogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgCgptb250aGx5X2JheWVzX2ZpdCA9IGxhcHBseShtb250aGx5X21heCwgYmF5ZXNfZml0dGVyX3BhcmFtX3NlYXJjaCwgZG9fZGlhZ24gPSBGQUxTRSwgCiAgICAgICAgICAgICAgICAgICAgICAgICAgIGRvX2F1dG9yZWcgPSBGQUxTRSwKICAgICAgICAgICAgICAgICAgICAgICAgICAgcHNkID0gYyg1MDAsMC4zLDAuMyksIG5pdF9mdWxsPTMwMDAwLCBuaXRfc2VhcmNoID0gMzAwMDAsCiAgICAgICAgICAgICAgICAgICAgICAgICAgIHRoaW5uaW5nID0gMzAwKQphY2NlcHRhbmNlX3JhdGVzID0gbGFwcGx5KG1vbnRobHlfYmF5ZXNfZml0LCBmdW5jdGlvbih4KSB4JGFjY2VwdGFuY2VfcmF0ZVsxLCBdKQpwcmludChhY2NlcHRhbmNlX3JhdGVzKQpiYXllc19wYXJhbXMgPSBsYXBwbHkobW9udGhseV9iYXllc19maXQsIGZ1bmN0aW9uKHgpIGFwcGx5KHgkcG9zdGVyaW9yLCAyLCBtZWFuKSkKYmF5ZXNfc3RkZXJyID0gbGFwcGx5KG1vbnRobHlfYmF5ZXNfZml0LCBmdW5jdGlvbih4KSBhcHBseSh4JHBvc3RlcmlvciwgMiwgc2QpKQoKCmBgYAoKClRPRE8tPiBSIGxhcmdlc3QgZml0CgoKCioqUEFSVCAyKioKRmlyc3QsIHdlIGNoZWNrIGlmIHRoZSBsb2NhdGlvbiBwYXJhbWV0ZXIgZGVwZW5kcyBvbiB0aW1lIHVzaW5nIGEgbGlrZWxpaG9vZCByYXRpbyB0ZXN0CmBgYHtyfQojbW9udGhseV9maXRzID0gbGFwcGx5KG1vbnRobHlfbWF4LCAKIyAgICAgICAgICAgICAgICAgICAgICBmdW5jdGlvbih4KSBmZ2V2KGRhdGEuZnJhbWUoeClbLDFdLCBtZXRob2Q9Ik5lbGRlci1NZWFkIikpCnJhdGlvcyA9IGxpc3QoKQp0cmVuZCA9IDE6bGVuZ3RoKGFzLmRhdGEuZnJhbWUobW9udGhseV9tYXhbWzEyXV0pJFYxKQp0cmVuZCA9ICh0cmVuZC1tZWFuKHRyZW5kKSkvc2QodHJlbmQpICMgc2NhbGUgYW5kIGNlbnRlciBjb3ZhcmlhdGVzIGFzIHJlY29tbWVuZGVkCmVycm9yX2Nhc2VzID0gYyg5LCAxMikKZm9yIChpIGluIDE6bGVuZ3RoKG1vbnRobHlfbWF4KSkgewogIHByaW50KGkpCiAgCiAgZml0X2NvbnN0ID0gZmdldihhcy5kYXRhLmZyYW1lKG1vbnRobHlfbWF4W1tpXV0pJFYxLCAKICAgICAgICAgICAgICAgICAgICAgICAgICAgICBtZXRob2QgPSAiTmVsZGVyLU1lYWQiLAogICAgICAgICAgICAgICAgICAgICAgICAgICAgIHN0ZC5lcnIgPSBGQUxTRSkKICBmaXRfZGVwZW5kYW50ID0gZmdldihhcy5kYXRhLmZyYW1lKG1vbnRobHlfbWF4W1tpXV0pJFYxLCAKICAgICAgICAgICAgICAgICAgICAgICAgICAgICBtZXRob2QgPSAiTmVsZGVyLU1lYWQiLAogICAgICAgICAgICAgICAgICAgICAgICAgICAgIG5zbG9jID0gdHJlbmQsCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgc3RkLmVyciA9IEZBTFNFKQogIAogIHJhdGlvc1tbaV1dID0gZml0X2NvbnN0JGRldi1maXRfZGVwZW5kYW50JGRldiAKfQpuYW1lcyhyYXRpb3MpID0gbmFtZXMobW9udGhseV9tYXgpCmNoaV85NWxldmVsID0gcWNoaXNxKDEtMC4wNS8xMiwxKQoKcGxvdCh1bmxpc3QocmF0aW9zKSxtYWluPSI5NSUgY29uZmlkZW5jZSB0ZXN0IGZvciB0aW1lIGluZGVwZW5kYW5jZSwgXG5Cb25mZXJyb25pIG11bHRpcGxlIFRlc3RpbmciLCB4bGFiPSJNb250aCIseWxhYj0iTGlrZWxpaG9vZCBSYXRpbyBTdGF0aXN0aWMiKQphYmxpbmUoYT1jaGlfOTVsZXZlbCxiPTAsY29sPSJyZWQiKQoKYGBgCgpOb3csIGxldCdzIGNoZWNrIGZvciBpbmRlcGVuZGFuY2UgZnJvbSBFTlNPCmBgYHtyfQojbW9udGhseV9maXRzID0gbGFwcGx5KG1vbnRobHlfbWF4LCAKIyAgICAgICAgICAgICAgICAgICAgICBmdW5jdGlvbih4KSBmZ2V2KGRhdGEuZnJhbWUoeClbLDFdLCBtZXRob2Q9Ik5lbGRlci1NZWFkIikpCnJhdGlvcyA9IGxpc3QoKQojIHNwbGl0IG5pbm8gZGF0YSBpbnRvIG1vbnRocwpuID0gbmlubzM0CmRpbShuKT1jKDEyLGxlbmd0aChhcy5kYXRhLmZyYW1lKG1vbnRobHlfbWF4W1sxMl1dKSRWMSkpCmVycm9yX2Nhc2VzID0gYyg5LCAxMikKZm9yIChpIGluIDE6bGVuZ3RoKG1vbnRobHlfbWF4KSkgewogIHByaW50KGkpCiAgdHJlbmQgPSBuW2ksXQogIHRyZW5kID0gKHRyZW5kLW1lYW4odHJlbmQpKS9zZCh0cmVuZCkgIyBzY2FsZSBhbmQgY2VudGVyIGNvdmFyaWF0ZXMgYXMgcmVjb21tZW5kZWQKICBmaXRfY29uc3QgPSBmZ2V2KGFzLmRhdGEuZnJhbWUobW9udGhseV9tYXhbW2ldXSkkVjEsIAogICAgICAgICAgICAgICAgICAgICAgICAgICAgIG1ldGhvZCA9ICJOZWxkZXItTWVhZCIsCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgc3RkLmVyciA9IEZBTFNFKQogIGZpdF9kZXBlbmRhbnQgPSBmZ2V2KGFzLmRhdGEuZnJhbWUobW9udGhseV9tYXhbW2ldXSkkVjEsIAogICAgICAgICAgICAgICAgICAgICAgICAgICAgIG1ldGhvZCA9ICJOZWxkZXItTWVhZCIsCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgbnNsb2MgPSB0cmVuZCwKICAgICAgICAgICAgICAgICAgICAgICAgICAgICBzdGQuZXJyID0gRkFMU0UpCiAgCiAgcmF0aW9zW1tpXV0gPSBmaXRfY29uc3QkZGV2LWZpdF9kZXBlbmRhbnQkZGV2IAp9Cm5hbWVzKHJhdGlvcykgPSBuYW1lcyhtb250aGx5X21heCkKY2hpXzk1bGV2ZWwgPSBxY2hpc3EoMS0wLjA1LzEyLDEpCgpwbG90KHVubGlzdChyYXRpb3MpLG1haW49Ijk1JSBjb25maWRlbmNlIHRlc3QgZm9yIGluZGVwZW5kYW5jZSBmcm9tIEVOU08sIFxuQm9uZmVycm9uaSBtdWx0aXBsZSB0ZXN0aW5nIiwgeGxhYj0iTW9udGgiLHlsYWI9Ikxpa2VseWhvb2QgUmF0aW8gU3RhdGlzdGljIikKYWJsaW5lKGE9Y2hpXzk1bGV2ZWwsYj0wLGNvbD0icmVkIikKCmBgYAoKQW5vdGhlciBtZXRob2QgaXMgdGhlIGNoaSBwbG90OgpgYGB7cn0KIyBwbG90IHRoZSBjaGkgcGxvdCBmb3IgZGVwZW5kYW5jZSBhbmFseXNpcwp0cmVuZCA9IDE6bGVuZ3RoKGFzLmRhdGEuZnJhbWUobW9udGhseV9tYXhbWzEyXV0pJFYxKQpuID0gbmlubzM0CmRpbShuKT1jKDEyLGxlbmd0aChhcy5kYXRhLmZyYW1lKG1vbnRobHlfbWF4W1sxMl1dKSRWMSkpCgpmb3IgKGkgaW4gMTpsZW5ndGgobW9udGhseV9tYXgpKSB7CiAgbmlubyA9IG5baSxdCiAgbV9kYXRhPWFzLmRhdGEuZnJhbWUobW9udGhseV9tYXhbW2ldXSkkVjEKICBkYXQubTFfbW9udGggPSBjYmluZChtX2RhdGEsdHJlbmQpOwogIGRhdC5tMV9uaW5vID0gY2JpbmQobV9kYXRhLG5pbm8pOwogIHBhcihtZnJvdz1jKDIsMikpCiAgY2hpcGxvdChkYXQubTFfbW9udGgsbWFpbjEgPSAiQ2hpIFBsb3QgVGltZSIsbWFpbjIgPSAiQ2hpIEJhciBQbG90IFRpbWUiKTsKICBjaGlwbG90KGRhdC5tMV9uaW5vLG1haW4xID0gIkNoaSBQbG90IEVOU08iLG1haW4yID0gIkNoaSBCYXIgUGxvdCBFTlNPIik7Cn0KYGBgCgoKCioqUEFSVCAzKioKV2Ugd2lsbCBub3cgYW5hbHlzZSB0ZW1wb3JhbCBjbHVzdGVyaW5nIG9mIGV4dHJlbWVzLiBGb3IgdGhpcywgd2Ugd2lsbCB1c2UgdGhlIGV4aXBsb3QgZnVuY3Rpb24gZnJvbSB0aGUgZXZkIGxpYnJhcnkuCgpgYGB7cn0KIyBkZWZpbmUgYSBmdW5jdGlvbiBmb3IgZ2V0dGluZyB0aGUgZXh0cmVtYWwgaW5kaWNlcyBmb3IgZWFjaCBtb250aCBmb3IgYSBnaXZlbiB0aHJlc2hvbGQKbW9udGhseV9laW5kZXggPC0gZnVuY3Rpb24oZGF0YSwgdGhyZXNob2xkX3AsIHI9MCl7CiAgZWkgPSBsaXN0KCkKICBmb3IgKGkgaW4gMTpsZW5ndGgoZGF0YSkpIHsKICAgIHRocmVzaG9sZCA9IHF1YW50aWxlKGFzLmRhdGEuZnJhbWUoZGF0YVtbaV1dKSRWMSwgdGhyZXNob2xkX3ApCiAgICBlaVtbaV1dPWV4aShhcy5kYXRhLmZyYW1lKGRhdGFbW2ldXSkkVjEsIHRocmVzaG9sZCwgcikKICB9CiAgbmFtZXMoZWkpID0gbmFtZXMoZGF0YSkKCiAgcmV0dXJuKGVpKQp9CgplaSA9IG1vbnRobHlfZWluZGV4KGdldF9tb250aGx5KHByb2RfdHMpLCAwLjk1KQpwbG90KHVubGlzdChlaSksIG1haW49IkV4dHJlbWFsIEluZGV4IGJ5IE1vbnRoLCA5NSUtUXVhbnRpbGUgYXMgVGhyZXNob2xkIiwgeGxhYj0iTW9udGgiLCB5bGFiPSJFeHRyZW1hbCBpbmRleCIpCmBgYAoKV2Ugb2JzZXJ2ZSB0aGF0IHRoZSBleHRyZW1hbCBpbmRleCBpcyB+MC4yNS1+MC40NSwgd2UgY2FuIHRoZXJlZm9yZSBjb25jbHVkZSB0aGF0IHdlIGhhdmUgc3Ryb25nIGRlcGVuZGFuY2Ugb2YgZXh0cmVtZXMsIHdpdGggdGhlIGxpbWl0aW5nIG1lYW4gY2x1c3RlciBzaXplIGJlaW5nIHJvdWdobHkgZnJvbSAyIHRvIDQuIFRoZSBjbHVzdGVyaW5nIGhhcyBubyBlZmZlY3QgZm9yIGVzdGltYXRlcnMgYmFzZWQgb24gdGhlIChtb250aGx5KSBtYXhpbXVtLCBidXQgdGhlIHIgbGFyZ2VzdCBlc3RpbWF0b3IgaXMgaW5mbHVlbmNlZCBieSBpdC4KCgoqKlBBUlQgNCoqCkZpcnN0LCBsZXQncyBlc3RpbWF0ZSB0aGUgMTAgeWVhciByZXR1cm4gbGV2ZWwgdXNpbmcgcG9pbnQgcHJvY2VzcyBhcHByb2FjaApgYGB7cn0KbW9udGhseV9maXRzX3BwID0gbGlzdCgpCm1vbnRobHlfZGF0YSA9IGdldF9tb250aGx5KHByb2RfdHMpCmVycm9yX2Nhc2VzID0gYyg1LCA5LCApCm1vbnRoX2RheXMgPSBjKDMxLDI4LDMxLDMwLDMxLDMwLDMxLDMxLDMwLDMxLDMwLDMxKQpmb3IgKGkgaW4gMTpsZW5ndGgobW9udGhseV9tYXgpKSB7CiAgcHJpbnQoaSkKICB0aHJlc2hvbGQgPSBxdWFudGlsZShhcy5kYXRhLmZyYW1lKG1vbnRobHlfZGF0YVtbaV1dKSRWMSwgMC45NSkKICAKICBpZiAoaSAlaW4lIGVycm9yX2Nhc2VzKSB7CiAgICBtb250aGx5X2ZpdHNfcHBbW2ldXSA9IGZwb3QoYXMuZGF0YS5mcmFtZShtb250aGx5X2RhdGFbW2ldXSkkVjEsCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgdGhyZXNob2xkID0gdGhyZXNob2xkLAogICAgICAgICAgICAgICAgICAgICAgICAgICAgIG1vZGVsPSJwcCIsCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgbnBwID0gbW9udGhfZGF5c1tpXSo4LAogICAgICAgICAgICAgICAgICAgICAgICAgICAgIGNtYXggPSBUUlVFLAogICAgICAgICAgICAgICAgICAgICAgICAgICAgIHIgPSAxLAogICAgICAgICAgICAgICAgICAgICAgICAgICAgIHN0ZC5lcnIgPSBGQUxTRSwKICAgICAgICAgICAgICAgICAgICAgICAgICAgICBtZXRob2QgPSAiTmVsZGVyLU1lYWQiKQogIH0KICBlbHNlIHsKICAgIG1vbnRobHlfZml0c19wcFtbaV1dID0gZnBvdChhcy5kYXRhLmZyYW1lKG1vbnRobHlfZGF0YVtbaV1dKSRWMSwKICAgICAgICAgICAgICAgICAgICAgICAgICAgICB0aHJlc2hvbGQgPSB0aHJlc2hvbGQsCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgbW9kZWw9InBwIiwKICAgICAgICAgICAgICAgICAgICAgICAgICAgICBucHAgPSBtb250aF9kYXlzW2ldKjgsCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgY21heCA9IFRSVUUsCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgciA9IDEsCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgbWV0aG9kID0gIk5lbGRlci1NZWFkIikKICB9Cn0KbmFtZXMobW9udGhseV9maXRzX3BwKSA9IG5hbWVzKG1vbnRobHlfZGF0YSkKZm9yKGkgaW4gMToxMil7CiAgcGFyKG1mcm93PWMoMiwyKSkgCiAgcGxvdChtb250aGx5X2ZpdHNfcHBbW2ldXSkKfQoKYGBgClRoZSBmaXQgaW4gZmVicnVhcnkgaGFzIGNvbXBsZXRlbHkgZmFpbGVkLCBhbmQgdGhlIG90aGVycyBhcmUgbm90IHZlcnkgZ29vZCBlaXRoZXIKCldlIHdpbGwgc3RpbGwgZXN0aW1hdGUgdGhlIHJldHVybiBsZXZlbHM6CmBgYHtyfQpyZXR1cm5fbGV2ZWwgPSBmdW5jdGlvbih4LHBlcmlvZD0yMCl7CiAgcCA9IDEvcGVyaW9kCiAgbG9jID0geCRlc3RpbWF0ZVtbMV1dCiAgc2NhbGUgPSB4JGVzdGltYXRlW1syXV0KICBzaGFwZSA9IHgkZXN0aW1hdGVbWzNdXQogIGxldmVsID0gbG9jICsgc2NhbGUqKCgoLWxvZygxLXApKV4tc2hhcGUtMSkvc2hhcGUpCiAgcmV0dXJuKGxldmVsKQp9CnJldHVybl9sZXZlbF8yMCA9IGxhcHBseShtb250aGx5X2ZpdHMsIHJldHVybl9sZXZlbCkgIyAyMCBmb3IgdGVzdGluZwpyZXR1cm5fbGV2ZWxfMTAwID0gbGFwcGx5KG1vbnRobHlfZml0cywgcmV0dXJuX2xldmVsLCBwZXJpb2Q9MTAwKQpyZXR1cm5fbGV2ZWxfNTAgPSBsYXBwbHkobW9udGhseV9maXRzLCByZXR1cm5fbGV2ZWwsIHBlcmlvZD01MCkKcGxvdCh1bmxpc3QocmV0dXJuX2xldmVsXzEwMCksbWFpbj0iMTAwIFllYXIgUmV0dXJuIGxldmVsLCBlc3RpbWF0ZWQgd2l0aCBwb2ludCBwcm9jZXNzIiwgeGxhYj0iTW9udGgiLHlsYWI9IlJldHVybiBMZXZlbCIpCnBsb3QodW5saXN0KHJldHVybl9sZXZlbF81MCksbWFpbj0iNTAgWWVhciBSZXR1cm4gbGV2ZWwsIGVzdGltYXRlZCB3aXRoIHBvaW50IHByb2Nlc3MiLCB4bGFiPSJNb250aCIseWxhYj0iUmV0dXJuIExldmVsIikKCmBgYAoKClRPRE86IGVzdGltYXQgd2l0aCBtY21jCkFzc3VtaW5nIHRoYXQgd2UgaGF2ZSB0aGUgcG9zdGVyaW9yIGRlbnNpdGllcyBvZiB0aGUgbWFya292IGNoYWlucywgY2FsbCB0aGVpcyBmdW5jdGlvbiB0byBwbG90IHJldHVybiBsZXZlbCBoaXN0b2dyYW1zCgpgYGB7cn0KcmV0dXJuX2xldmVsX21jbWMgPSBmdW5jdGlvbihwb3N0ZXJpb3IscGVyaW9kPTIwKXsKICB1ID0gbWMucXVhbnQocG9zdGVyaW9yLHA9MS0xL3BlcmlvZCxsaD0iZ2V2IikKICBsYWJlbF9tY21jX3JsID0gc3ByaW50ZigiJXMgWWVhciBSZXR1cm4gTGV2ZWwiLHBlcmlvZCkKICBoaXN0KHUsbmNsYXNzPTIwLHByb2I9VCx4bGFiPWxhYmVsX21jbWNfcmwsIG1haW4gPSAiUmV0dXJuIExldmVsIEhpc3RvZ3JhbSIpCn0KCmxhcHBseShtb250aGx5X2JheWVzX2ZpdCwgZnVuY3Rpb24oeCkgcmV0dXJuX2xldmVsX21jbWMoeCRwb3N0ZXJpb3IpKQoKCmBgYAoKCgoKCgoKCg==